Load packages

library(tidyverse)
library(gridExtra)
library(cowplot)
library(viridis)
library(ggridges)
library(ggstance)
library(treeio)
library(ggtree)
library(tidytree)

Load and combine data

studies <- read_csv("../data/studies_gsheet.csv")
species <- read_csv("../data/species_gsheet.csv")
fulltree <- read.nexus("../data/consensusTree_10kTrees_298Primates_V3.nex")
refs <- read_csv("../data/ref_nodes.csv")
data <- studies %>% 
  left_join(species, by = c("species" = "species_file")) %>% 
  rename(label = species_tree, num = n) %>% 
  select(label, studyID, species, site, num)
data2 <- data %>% group_by(label, species) %>% summarise(totalN = sum(num, na.rm = T))
# turn tree into tidy dataframe
tree2 <- fulltree %>% 
  drop.tip(c("Pan_troglodytes_schweinfurthii", "Pan_troglodytes_troglodytes",
             "Pan_troglodytes_vellerosus", "Pongo_pygmaeus")) %>%
  as_tibble

inner_nodes <- c(295:304, 314:315, 320:326, 328:332, 334, 336:353, 356:358, 362:364, 366, 382:390, 395, 399:402, 404:410, 412:418, 420:423, 430:433, 435:439, 445:449, 453:455, 457:458, 459:461, 464:471, 473:478, 480:484, 488:499, 512, 531:532, 536:538, 541:545, 552:554, 561:563, 566, 568)

tree3 <- tree2 %>% 
  mutate(label = fct_recode(label, "Pan_troglodytes" = "Pan_troglodytes_verus",
                                   "Pongo_spp" = "Pongo_abelii")) %>% 
  left_join(data2) %>% 
  mutate(
    hasN = ifelse(is.na(totalN), .1, .5), # used to size branches + color the tip labels
    hasN2 = ifelse(is.na(totalN) & !(node %in% inner_nodes), 0, .5), # used to color branches
    label = str_replace_all(label, "_", " ")) %>% 
  left_join(refs) %>% 
  groupClade(refs$node[-1]) %>% 
  mutate(group = fct_recode(group, "2" = "1"))
Joining, by = "label"
Column `label` joining factor and character vector, coercing into character vectorJoining, by = "node"
# turn back into tree
tree4 <- as.treedata(tree3)

Circular tree of the 10ktree primates

cols <- viridis(4, end = .9)
p <- ggtree(tree4, aes(size = hasN, alpha = hasN2), layout = "circular") +
  # highlight clades with background colors
  geom_hilight(node = 485, fill = cols[1], alpha = .3) +
  geom_hilight(node = 488, fill = cols[1], alpha = .3) +
  geom_hilight(node = 421, fill = cols[2], alpha = .3) +
  geom_hilight(node = 299, fill = cols[3], alpha = .3) +
  geom_hilight(node = 404, fill = cols[4], alpha = .3) +
  # plot tree again to be on top of the highlights
  geom_tree() +
  # root
  geom_rootpoint(size = 1) +
  # tips
  geom_tippoint(aes(size = totalN), alpha = .7) +
  geom_tiplab2(aes(alpha = hasN), offset = 3, size = 3) +
  # tweak scales
  scale_alpha_continuous(range = c(.2, 1)) +
  scale_size_continuous(range = c(.5, 15)) +
  # widen plotting area
  xlim(NA, 100)
p

ggsave("../graphs/phylo_full.pdf", p, width = 8, height = 8, scale = 2)
# to figure out node numbers
n1 <- p + geom_text(aes(label = node, x = branch), size = 2, col = "blue", vjust = -.5)
ggsave("../graphs/full_tree_nodes_circular.pdf", n1, width = 8, height = 8, scale = 2)
n2 <- ggtree(tree4, aes(size = hasN, alpha = hasN2)) +
  # highlight clades with background colors
  geom_hilight(node = 485, fill = cols[1], alpha = .3) +
  geom_hilight(node = 488, fill = cols[1], alpha = .3) +
  geom_hilight(node = 421, fill = cols[2], alpha = .3) +
  geom_hilight(node = 299, fill = cols[3], alpha = .3) +
  geom_hilight(node = 404, fill = cols[4], alpha = .3) +
  # plot tree again to be on top of the highlights
  geom_tree() +
  # root
  geom_rootedge(rootedge = 2) +
  geom_rootpoint(size = 1) +
  # node labels
  geom_text(aes(label = node, x = branch), size = 2, col = "blue", vjust = -.5) +
  # tips
  geom_tippoint(aes(size = totalN), alpha = .7) +
  geom_tiplab(aes(alpha = hasN), offset = 1.8, size = 3) +
  # tweak scales
  scale_alpha_continuous(range = c(.2, 1)) +
  scale_size_continuous(range = c(.5, 15)) +
  # widen plotting area
  expand_limits(x = 90) +
  theme_tree2()

ggsave("../graphs/full_tree_nodes.pdf", n2, width = 8, height = 20, scale = 2)

Sample size in detail

data %>% arrange(desc(num)) %>% mutate(site = str_sub(site, 1, 35)) %>% select(-species)
# subset tree to just those species who have sample sizes reported, i.e. those who were tested
to_drop <- tree3 %>% filter(is.na(totalN)) %>% pull(label)
tree5 <- drop.tip(tree4, to_drop)
d3 <- mutate(data, label = str_replace_all(label, "_", " "))
# filter super large samples out for visualization? note in caption
# species with more than X sites can get a density
d3a <- d3 %>% group_by(species) %>% filter(n_distinct(site) >= 4, num <= 200)
d3b <- d3 %>% # setdiff(d3, d3a) %>% ## <- to NOT show points for densities
  group_by(species) %>% 
  # create variable num2 is NA if there's only one data point for a species
  # --> those species will only get the vertical crossbar
  mutate(flag = n_distinct(site) == 1) %>% 
  ungroup %>% 
  mutate(num2 = ifelse(flag, NA, num)) %>% 
  filter(num <= 200)

# for vertical crossbar = median
d4a <- d3 %>% 
  group_by(label, species) %>% 
  summarise(Mdn = median(num, na.rm = T)) # totalN = sum(num), sitesN = n_distinct(site)

# for cross = median (sample size by species and site)
d4b <- d3 %>% 
  group_by(label, species, site) %>% 
  summarise(Mdn_site = median(num, na.rm = T))%>% 
  group_by(label, species) %>% 
  summarise(Mdn_site2 = median(Mdn_site))

# make NA when closer than 1 to Mdn
d4 <- full_join(d4a, d4b) %>% 
  mutate(Mdn_site2 = ifelse(abs(Mdn - Mdn_site2) < 1, NA, Mdn_site2))
Joining, by = c("label", "species")
# for vertical line in ridge plot (grand median)
# + hacky way to make horizontal grid lines for right panel only
v <- tibble(reference = c(NA, median(d3$num, na.rm = T)), .panel = c("Tree", "xSample size"))
h <- tibble(reference = c(NA, 1:Ntip(tree5)), .panel = c("Tree", rep("xSample size", Ntip(tree5))))

# for axis labels
ax <- tibble(lab = c("Distance (Millions of years)", "Sample size"), 
             x = c(60, 100), y = -4, .panel = c("Tree", "xSample size"))
# LEFT FACET
q <- ggtree(tree5, aes(col = group)) +
  # root
  geom_rootedge(rootedge = 5) +
  # tip labels
  geom_tippoint(aes(size = totalN), alpha = .5) +
  geom_tiplab(offset = 4, size = 3) +
  # tweak scales
  scale_color_manual(values = c("grey30", cols)) +
  scale_fill_manual(values = cols) + # when all categories are taken: cols
  # display timescale at the bottom
  theme_tree2() +
  xlim_tree(112) +
  xlim_expand(c(0, 175), "xSample size") +
  # add axis labels
  geom_text(data = ax, aes(label = lab), col = "black") +
  scale_x_continuous(expand = expand_scale(mult = c(0, .01))) +
  scale_y_continuous(limits = c(2, Ntip(tree5)), oob = function(x, ...) x) +
  coord_cartesian(clip = "off") +
  # add reference lines (these will show up on right panel of facet_plot only)
  geom_hline(data = h, aes(yintercept = reference), lwd = .2, col = "grey", alpha = .5) +
  geom_vline(data = v, aes(xintercept = reference), lwd = 1.5, col = "grey", alpha = .3) +
  # remove facet strips, expand bottom margin (to make space for x axis labels)
  theme(strip.text = element_blank(), strip.background = element_blank(),
        plot.margin = unit(c(1, 1, 2, 1.5), "cm"), panel.spacing = unit(1, "cm"))
# right-side viz depends on the number of sites per species:
# 1 site = vertical crossbar only
# 2+ sites = points + crossbar at median
# X+ sites = densities (currently, X = 4 just to illustrate)

# dirty hack: x in front of "Sample size" is to have that panel sort to the right (alphabetically) until I figure out why it doesn't just go by order. This cropped up as an issue when I added the dummy point for the x-axis expansion...

# ADD RIGHT FACET
q %>% 
  # densities for species with enough sites
  facet_plot("xSample size", d3a, geom_density_ridges, 
             aes(x = num, group = label, fill = group, height = ..density..),
             alpha = .5, lwd = .3, scale = .3) %>%
  # vertical crossbar for Mdn, X for Mdn of site medians
  facet_plot("xSample size", d4, geom_crossbarh, aes(x = Mdn, xmin = Mdn, xmax = Mdn, group = label,
             col = group), alpha = .5, width = .6, fatten = 1.5) %>%
  facet_plot("xSample size", d4, geom_point, aes(x = Mdn_site2, group = label), shape = 4, size = 2.5,
             stroke = 1.05) %>%
  # vertical mark for individual sites
  facet_plot("xSample size", d3b, geom_jitter, aes(x = num2, group = label), shape = "|", size = 2.5,
             width = .1, height = 0, alpha = .8)

ggsave("../graphs/phylo_ridge_site.pdf", width = 6, height = 6, scale = 2)
# ADD RIGHT FACET
q %>% 
  facet_plot("xSample size", filter(d3, num < 200), geom_boxploth, aes(x = num, group = label, 
             fill = group), alpha = .5, lwd = .3, varwidth = T) %>% 
  # vertical crossbar for Mdn
  facet_plot("xSample size", d4, geom_crossbarh, aes(x = Mdn, xmin = Mdn, xmax = Mdn, group = label,
             col = group), alpha = .5, width = .6, fatten = 1.5)

ggsave("../graphs/phylo_ridge_box.pdf", width = 6, height = 8, scale = 2)

Session info

sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Mojave 10.14.5

Matrix products: default
BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] tidytree_0.2.8    ggtree_1.16.6     treeio_1.8.2      ggstance_0.3.3    ggridges_0.5.1   
 [6] viridis_0.5.1     viridisLite_0.3.0 cowplot_1.0.0     gridExtra_2.3     forcats_0.4.0    
[11] stringr_1.4.0     dplyr_0.8.3       purrr_0.3.2       readr_1.3.1       tidyr_1.0.0      
[16] tibble_2.1.3      ggplot2_3.2.1     tidyverse_1.2.1  

loaded via a namespace (and not attached):
 [1] tidyselect_0.2.5   xfun_0.10          reshape2_1.4.3     haven_2.1.1        lattice_0.20-38   
 [6] colorspace_1.4-1   vctrs_0.2.0        generics_0.0.2     htmltools_0.4.0    base64enc_0.1-3   
[11] yaml_2.2.0         rlang_0.4.0        pillar_1.4.2       glue_1.3.1         withr_2.1.2       
[16] modelr_0.1.5       readxl_1.3.1       rvcheck_0.1.5      lifecycle_0.1.0    plyr_1.8.4        
[21] munsell_0.5.0      gtable_0.3.0       cellranger_1.1.0   rvest_0.3.4        evaluate_0.14     
[26] labeling_0.3       knitr_1.25         parallel_3.6.1     broom_0.5.2        Rcpp_1.0.2        
[31] BiocManager_1.30.7 scales_1.0.0       backports_1.1.5    jsonlite_1.6       digest_0.6.21     
[36] hms_0.5.1          stringi_1.4.3      grid_3.6.1         cli_1.1.0          tools_3.6.1       
[41] magrittr_1.5       lazyeval_0.2.2     crayon_1.3.4       ape_5.3            pkgconfig_2.0.3   
[46] zeallot_0.1.0      xml2_1.2.2         lubridate_1.7.4    rmarkdown_1.16     assertthat_0.2.1  
[51] httr_1.4.1         rstudioapi_0.10    R6_2.4.0           nlme_3.1-141       compiler_3.6.1    
LS0tCnRpdGxlOiAiUGh5bG9nZW5ldGljIFRyZWUiCm91dHB1dDogCiAgaHRtbF9ub3RlYm9vazoKICAgIGNzczogc3R5bGUuY3NzCiAgICB0aGVtZTogcGFwZXIKLS0tCgpMb2FkIHBhY2thZ2VzCgpgYGB7ciwgbWVzc2FnZT1GQUxTRX0KbGlicmFyeSh0aWR5dmVyc2UpCmxpYnJhcnkoZ3JpZEV4dHJhKQpsaWJyYXJ5KGNvd3Bsb3QpCmxpYnJhcnkodmlyaWRpcykKbGlicmFyeShnZ3JpZGdlcykKbGlicmFyeShnZ3N0YW5jZSkKbGlicmFyeSh0cmVlaW8pCmxpYnJhcnkoZ2d0cmVlKQpsaWJyYXJ5KHRpZHl0cmVlKQpgYGAKCkxvYWQgYW5kIGNvbWJpbmUgZGF0YQoKYGBge3IsIG1lc3NhZ2U9RkFMU0V9CnN0dWRpZXMgPC0gcmVhZF9jc3YoIi4uL2RhdGEvc3R1ZGllc19nc2hlZXQuY3N2IikKc3BlY2llcyA8LSByZWFkX2NzdigiLi4vZGF0YS9zcGVjaWVzX2dzaGVldC5jc3YiKQpmdWxsdHJlZSA8LSByZWFkLm5leHVzKCIuLi9kYXRhL2NvbnNlbnN1c1RyZWVfMTBrVHJlZXNfMjk4UHJpbWF0ZXNfVjMubmV4IikKcmVmcyA8LSByZWFkX2NzdigiLi4vZGF0YS9yZWZfbm9kZXMuY3N2IikKYGBgCgpgYGB7cn0KZGF0YSA8LSBzdHVkaWVzICU+JSAKICBsZWZ0X2pvaW4oc3BlY2llcywgYnkgPSBjKCJzcGVjaWVzIiA9ICJzcGVjaWVzX2ZpbGUiKSkgJT4lIAogIHJlbmFtZShsYWJlbCA9IHNwZWNpZXNfdHJlZSwgbnVtID0gbikgJT4lIAogIHNlbGVjdChsYWJlbCwgc3R1ZHlJRCwgc3BlY2llcywgc2l0ZSwgbnVtKQpgYGAKCmBgYHtyfQpkYXRhMiA8LSBkYXRhICU+JSBncm91cF9ieShsYWJlbCwgc3BlY2llcykgJT4lIHN1bW1hcmlzZSh0b3RhbE4gPSBzdW0obnVtLCBuYS5ybSA9IFQpKQpgYGAKCmBgYHtyfQojIHR1cm4gdHJlZSBpbnRvIHRpZHkgZGF0YWZyYW1lCnRyZWUyIDwtIGZ1bGx0cmVlICU+JSAKICBkcm9wLnRpcChjKCJQYW5fdHJvZ2xvZHl0ZXNfc2Nod2VpbmZ1cnRoaWkiLCAiUGFuX3Ryb2dsb2R5dGVzX3Ryb2dsb2R5dGVzIiwKICAgICAgICAgICAgICJQYW5fdHJvZ2xvZHl0ZXNfdmVsbGVyb3N1cyIsICJQb25nb19weWdtYWV1cyIpKSAlPiUKICBhc190aWJibGUKCmlubmVyX25vZGVzIDwtIGMoMjk1OjMwNCwgMzE0OjMxNSwgMzIwOjMyNiwgMzI4OjMzMiwgMzM0LCAzMzY6MzUzLCAzNTY6MzU4LCAzNjI6MzY0LCAzNjYsIDM4MjozOTAsIDM5NSwgMzk5OjQwMiwgNDA0OjQxMCwgNDEyOjQxOCwgNDIwOjQyMywgNDMwOjQzMywgNDM1OjQzOSwgNDQ1OjQ0OSwgNDUzOjQ1NSwgNDU3OjQ1OCwgNDU5OjQ2MSwgNDY0OjQ3MSwgNDczOjQ3OCwgNDgwOjQ4NCwgNDg4OjQ5OSwgNTEyLCA1MzE6NTMyLCA1MzY6NTM4LCA1NDE6NTQ1LCA1NTI6NTU0LCA1NjE6NTYzLCA1NjYsIDU2OCkKCnRyZWUzIDwtIHRyZWUyICU+JSAKICBtdXRhdGUobGFiZWwgPSBmY3RfcmVjb2RlKGxhYmVsLCAiUGFuX3Ryb2dsb2R5dGVzIiA9ICJQYW5fdHJvZ2xvZHl0ZXNfdmVydXMiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJQb25nb19zcHAiID0gIlBvbmdvX2FiZWxpaSIpKSAlPiUgCiAgbGVmdF9qb2luKGRhdGEyKSAlPiUgCiAgbXV0YXRlKAogICAgaGFzTiA9IGlmZWxzZShpcy5uYSh0b3RhbE4pLCAuMSwgLjUpLCAjIHVzZWQgdG8gc2l6ZSBicmFuY2hlcyArIGNvbG9yIHRoZSB0aXAgbGFiZWxzCiAgICBoYXNOMiA9IGlmZWxzZShpcy5uYSh0b3RhbE4pICYgIShub2RlICVpbiUgaW5uZXJfbm9kZXMpLCAwLCAuNSksICMgdXNlZCB0byBjb2xvciBicmFuY2hlcwogICAgbGFiZWwgPSBzdHJfcmVwbGFjZV9hbGwobGFiZWwsICJfIiwgIiAiKSkgJT4lIAogIGxlZnRfam9pbihyZWZzKSAlPiUgCiAgZ3JvdXBDbGFkZShyZWZzJG5vZGVbLTFdKSAlPiUgCiAgbXV0YXRlKGdyb3VwID0gZmN0X3JlY29kZShncm91cCwgIjIiID0gIjEiKSkKCiMgdHVybiBiYWNrIGludG8gdHJlZQp0cmVlNCA8LSBhcy50cmVlZGF0YSh0cmVlMykKYGBgCgojIENpcmN1bGFyIHRyZWUgb2YgdGhlIDEwa3RyZWUgcHJpbWF0ZXMKCmBgYHtyfQpjb2xzIDwtIHZpcmlkaXMoNCwgZW5kID0gLjkpCmBgYAoKYGBge3J9CnAgPC0gZ2d0cmVlKHRyZWU0LCBhZXMoc2l6ZSA9IGhhc04sIGFscGhhID0gaGFzTjIpLCBsYXlvdXQgPSAiY2lyY3VsYXIiKSArCiAgIyBoaWdobGlnaHQgY2xhZGVzIHdpdGggYmFja2dyb3VuZCBjb2xvcnMKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQ4NSwgZmlsbCA9IGNvbHNbMV0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQ4OCwgZmlsbCA9IGNvbHNbMV0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQyMSwgZmlsbCA9IGNvbHNbMl0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDI5OSwgZmlsbCA9IGNvbHNbM10sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQwNCwgZmlsbCA9IGNvbHNbNF0sIGFscGhhID0gLjMpICsKICAjIHBsb3QgdHJlZSBhZ2FpbiB0byBiZSBvbiB0b3Agb2YgdGhlIGhpZ2hsaWdodHMKICBnZW9tX3RyZWUoKSArCiAgIyByb290CiAgZ2VvbV9yb290cG9pbnQoc2l6ZSA9IDEpICsKICAjIHRpcHMKICBnZW9tX3RpcHBvaW50KGFlcyhzaXplID0gdG90YWxOKSwgYWxwaGEgPSAuNykgKwogIGdlb21fdGlwbGFiMihhZXMoYWxwaGEgPSBoYXNOKSwgb2Zmc2V0ID0gMywgc2l6ZSA9IDMpICsKICAjIHR3ZWFrIHNjYWxlcwogIHNjYWxlX2FscGhhX2NvbnRpbnVvdXMocmFuZ2UgPSBjKC4yLCAxKSkgKwogIHNjYWxlX3NpemVfY29udGludW91cyhyYW5nZSA9IGMoLjUsIDE1KSkgKwogICMgd2lkZW4gcGxvdHRpbmcgYXJlYQogIHhsaW0oTkEsIDEwMCkKYGBgCgpgYGB7ciwgZmlnLndpZHRoPTgsIGZpZy5oZWlnaHQ9OCwgY2FjaGU9VFJVRX0KcApgYGAKCmBgYHtyLCBjYWNoZT1UUlVFfQpnZ3NhdmUoIi4uL2dyYXBocy9waHlsb19mdWxsLnBkZiIsIHAsIHdpZHRoID0gOCwgaGVpZ2h0ID0gOCwgc2NhbGUgPSAyKQpgYGAKCmBgYHtyLCBmaWcud2lkdGg9OCwgZmlnLmhlaWdodD04LCBjYWNoZT1UUlVFLCBldmFsPUZBTFNFfQojIHRvIGZpZ3VyZSBvdXQgbm9kZSBudW1iZXJzCm4xIDwtIHAgKyBnZW9tX3RleHQoYWVzKGxhYmVsID0gbm9kZSwgeCA9IGJyYW5jaCksIHNpemUgPSAyLCBjb2wgPSAiYmx1ZSIsIHZqdXN0ID0gLS41KQpnZ3NhdmUoIi4uL2dyYXBocy9mdWxsX3RyZWVfbm9kZXNfY2lyY3VsYXIucGRmIiwgbjEsIHdpZHRoID0gOCwgaGVpZ2h0ID0gOCwgc2NhbGUgPSAyKQpgYGAKCmBgYHtyLCBmaWcud2lkdGg9OCwgZmlnLmhlaWdodD0yMCwgY2FjaGU9VFJVRSwgZXZhbD1GQUxTRX0KbjIgPC0gZ2d0cmVlKHRyZWU0LCBhZXMoc2l6ZSA9IGhhc04sIGFscGhhID0gaGFzTjIpKSArCiAgIyBoaWdobGlnaHQgY2xhZGVzIHdpdGggYmFja2dyb3VuZCBjb2xvcnMKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQ4NSwgZmlsbCA9IGNvbHNbMV0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQ4OCwgZmlsbCA9IGNvbHNbMV0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQyMSwgZmlsbCA9IGNvbHNbMl0sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDI5OSwgZmlsbCA9IGNvbHNbM10sIGFscGhhID0gLjMpICsKICBnZW9tX2hpbGlnaHQobm9kZSA9IDQwNCwgZmlsbCA9IGNvbHNbNF0sIGFscGhhID0gLjMpICsKICAjIHBsb3QgdHJlZSBhZ2FpbiB0byBiZSBvbiB0b3Agb2YgdGhlIGhpZ2hsaWdodHMKICBnZW9tX3RyZWUoKSArCiAgIyByb290CiAgZ2VvbV9yb290ZWRnZShyb290ZWRnZSA9IDIpICsKICBnZW9tX3Jvb3Rwb2ludChzaXplID0gMSkgKwogICMgbm9kZSBsYWJlbHMKICBnZW9tX3RleHQoYWVzKGxhYmVsID0gbm9kZSwgeCA9IGJyYW5jaCksIHNpemUgPSAyLCBjb2wgPSAiYmx1ZSIsIHZqdXN0ID0gLS41KSArCiAgIyB0aXBzCiAgZ2VvbV90aXBwb2ludChhZXMoc2l6ZSA9IHRvdGFsTiksIGFscGhhID0gLjcpICsKICBnZW9tX3RpcGxhYihhZXMoYWxwaGEgPSBoYXNOKSwgb2Zmc2V0ID0gMS44LCBzaXplID0gMykgKwogICMgdHdlYWsgc2NhbGVzCiAgc2NhbGVfYWxwaGFfY29udGludW91cyhyYW5nZSA9IGMoLjIsIDEpKSArCiAgc2NhbGVfc2l6ZV9jb250aW51b3VzKHJhbmdlID0gYyguNSwgMTUpKSArCiAgIyB3aWRlbiBwbG90dGluZyBhcmVhCiAgZXhwYW5kX2xpbWl0cyh4ID0gOTApICsKICB0aGVtZV90cmVlMigpCgpnZ3NhdmUoIi4uL2dyYXBocy9mdWxsX3RyZWVfbm9kZXMucGRmIiwgbjIsIHdpZHRoID0gOCwgaGVpZ2h0ID0gMjAsIHNjYWxlID0gMikKYGBgCgoKIyBTYW1wbGUgc2l6ZSBpbiBkZXRhaWwKCmBgYHtyfQpkYXRhICU+JSBhcnJhbmdlKGRlc2MobnVtKSkgJT4lIG11dGF0ZShzaXRlID0gc3RyX3N1YihzaXRlLCAxLCAzNSkpICU+JSBzZWxlY3QoLXNwZWNpZXMpCmBgYAoKYGBge3J9CiMgc3Vic2V0IHRyZWUgdG8ganVzdCB0aG9zZSBzcGVjaWVzIHdobyBoYXZlIHNhbXBsZSBzaXplcyByZXBvcnRlZCwgaS5lLiB0aG9zZSB3aG8gd2VyZSB0ZXN0ZWQKdG9fZHJvcCA8LSB0cmVlMyAlPiUgZmlsdGVyKGlzLm5hKHRvdGFsTikpICU+JSBwdWxsKGxhYmVsKQp0cmVlNSA8LSBkcm9wLnRpcCh0cmVlNCwgdG9fZHJvcCkKZDMgPC0gbXV0YXRlKGRhdGEsIGxhYmVsID0gc3RyX3JlcGxhY2VfYWxsKGxhYmVsLCAiXyIsICIgIikpCmBgYAoKYGBge3J9CiMgZmlsdGVyIHN1cGVyIGxhcmdlIHNhbXBsZXMgb3V0IGZvciB2aXN1YWxpemF0aW9uPyBub3RlIGluIGNhcHRpb24KIyBzcGVjaWVzIHdpdGggbW9yZSB0aGFuIFggc2l0ZXMgY2FuIGdldCBhIGRlbnNpdHkKZDNhIDwtIGQzICU+JSBncm91cF9ieShzcGVjaWVzKSAlPiUgZmlsdGVyKG5fZGlzdGluY3Qoc2l0ZSkgPj0gNCwgbnVtIDw9IDIwMCkKZDNiIDwtIGQzICU+JSAjIHNldGRpZmYoZDMsIGQzYSkgJT4lICMjIDwtIHRvIE5PVCBzaG93IHBvaW50cyBmb3IgZGVuc2l0aWVzCiAgZ3JvdXBfYnkoc3BlY2llcykgJT4lIAogICMgY3JlYXRlIHZhcmlhYmxlIG51bTIgaXMgTkEgaWYgdGhlcmUncyBvbmx5IG9uZSBkYXRhIHBvaW50IGZvciBhIHNwZWNpZXMKICAjIC0tPiB0aG9zZSBzcGVjaWVzIHdpbGwgb25seSBnZXQgdGhlIHZlcnRpY2FsIGNyb3NzYmFyCiAgbXV0YXRlKGZsYWcgPSBuX2Rpc3RpbmN0KHNpdGUpID09IDEpICU+JSAKICB1bmdyb3VwICU+JSAKICBtdXRhdGUobnVtMiA9IGlmZWxzZShmbGFnLCBOQSwgbnVtKSkgJT4lIAogIGZpbHRlcihudW0gPD0gMjAwKQoKIyBmb3IgdmVydGljYWwgY3Jvc3NiYXIgPSBtZWRpYW4KZDRhIDwtIGQzICU+JSAKICBncm91cF9ieShsYWJlbCwgc3BlY2llcykgJT4lIAogIHN1bW1hcmlzZShNZG4gPSBtZWRpYW4obnVtLCBuYS5ybSA9IFQpKSAjIHRvdGFsTiA9IHN1bShudW0pLCBzaXRlc04gPSBuX2Rpc3RpbmN0KHNpdGUpCgojIGZvciBjcm9zcyA9IG1lZGlhbiAoc2FtcGxlIHNpemUgYnkgc3BlY2llcyBhbmQgc2l0ZSkKZDRiIDwtIGQzICU+JSAKICBncm91cF9ieShsYWJlbCwgc3BlY2llcywgc2l0ZSkgJT4lIAogIHN1bW1hcmlzZShNZG5fc2l0ZSA9IG1lZGlhbihudW0sIG5hLnJtID0gVCkpJT4lIAogIGdyb3VwX2J5KGxhYmVsLCBzcGVjaWVzKSAlPiUgCiAgc3VtbWFyaXNlKE1kbl9zaXRlMiA9IG1lZGlhbihNZG5fc2l0ZSkpCgojIG1ha2UgTkEgd2hlbiBjbG9zZXIgdGhhbiAxIHRvIE1kbgpkNCA8LSBmdWxsX2pvaW4oZDRhLCBkNGIpICU+JSAKICBtdXRhdGUoTWRuX3NpdGUyID0gaWZlbHNlKGFicyhNZG4gLSBNZG5fc2l0ZTIpIDwgMSwgTkEsIE1kbl9zaXRlMikpCgojIGZvciB2ZXJ0aWNhbCBsaW5lIGluIHJpZGdlIHBsb3QgKGdyYW5kIG1lZGlhbikKIyArIGhhY2t5IHdheSB0byBtYWtlIGhvcml6b250YWwgZ3JpZCBsaW5lcyBmb3IgcmlnaHQgcGFuZWwgb25seQp2IDwtIHRpYmJsZShyZWZlcmVuY2UgPSBjKE5BLCBtZWRpYW4oZDMkbnVtLCBuYS5ybSA9IFQpKSwgLnBhbmVsID0gYygiVHJlZSIsICJ4U2FtcGxlIHNpemUiKSkKaCA8LSB0aWJibGUocmVmZXJlbmNlID0gYyhOQSwgMTpOdGlwKHRyZWU1KSksIC5wYW5lbCA9IGMoIlRyZWUiLCByZXAoInhTYW1wbGUgc2l6ZSIsIE50aXAodHJlZTUpKSkpCgojIGZvciBheGlzIGxhYmVscwpheCA8LSB0aWJibGUobGFiID0gYygiRGlzdGFuY2UgKE1pbGxpb25zIG9mIHllYXJzKSIsICJTYW1wbGUgc2l6ZSIpLCAKICAgICAgICAgICAgIHggPSBjKDYwLCAxMDApLCB5ID0gLTQsIC5wYW5lbCA9IGMoIlRyZWUiLCAieFNhbXBsZSBzaXplIikpCmBgYAoKYGBge3IsIGNhY2hlPVRSVUV9CiMgTEVGVCBGQUNFVApxIDwtIGdndHJlZSh0cmVlNSwgYWVzKGNvbCA9IGdyb3VwKSkgKwogICMgcm9vdAogIGdlb21fcm9vdGVkZ2Uocm9vdGVkZ2UgPSA1KSArCiAgIyB0aXAgbGFiZWxzCiAgZ2VvbV90aXBwb2ludChhZXMoc2l6ZSA9IHRvdGFsTiksIGFscGhhID0gLjUpICsKICBnZW9tX3RpcGxhYihvZmZzZXQgPSA0LCBzaXplID0gMykgKwogICMgdHdlYWsgc2NhbGVzCiAgc2NhbGVfY29sb3JfbWFudWFsKHZhbHVlcyA9IGMoImdyZXkzMCIsIGNvbHMpKSArCiAgc2NhbGVfZmlsbF9tYW51YWwodmFsdWVzID0gY29scykgKyAjIHdoZW4gYWxsIGNhdGVnb3JpZXMgYXJlIHRha2VuOiBjb2xzCiAgIyBkaXNwbGF5IHRpbWVzY2FsZSBhdCB0aGUgYm90dG9tCiAgdGhlbWVfdHJlZTIoKSArCiAgeGxpbV90cmVlKDExMikgKwogIHhsaW1fZXhwYW5kKGMoMCwgMTc1KSwgInhTYW1wbGUgc2l6ZSIpICsKICAjIGFkZCBheGlzIGxhYmVscwogIGdlb21fdGV4dChkYXRhID0gYXgsIGFlcyhsYWJlbCA9IGxhYiksIGNvbCA9ICJibGFjayIpICsKICBzY2FsZV94X2NvbnRpbnVvdXMoZXhwYW5kID0gZXhwYW5kX3NjYWxlKG11bHQgPSBjKDAsIC4wMSkpKSArCiAgc2NhbGVfeV9jb250aW51b3VzKGxpbWl0cyA9IGMoMiwgTnRpcCh0cmVlNSkpLCBvb2IgPSBmdW5jdGlvbih4LCAuLi4pIHgpICsKICBjb29yZF9jYXJ0ZXNpYW4oY2xpcCA9ICJvZmYiKSArCiAgIyBhZGQgcmVmZXJlbmNlIGxpbmVzICh0aGVzZSB3aWxsIHNob3cgdXAgb24gcmlnaHQgcGFuZWwgb2YgZmFjZXRfcGxvdCBvbmx5KQogIGdlb21faGxpbmUoZGF0YSA9IGgsIGFlcyh5aW50ZXJjZXB0ID0gcmVmZXJlbmNlKSwgbHdkID0gLjIsIGNvbCA9ICJncmV5IiwgYWxwaGEgPSAuNSkgKwogIGdlb21fdmxpbmUoZGF0YSA9IHYsIGFlcyh4aW50ZXJjZXB0ID0gcmVmZXJlbmNlKSwgbHdkID0gMS41LCBjb2wgPSAiZ3JleSIsIGFscGhhID0gLjMpICsKICAjIHJlbW92ZSBmYWNldCBzdHJpcHMsIGV4cGFuZCBib3R0b20gbWFyZ2luICh0byBtYWtlIHNwYWNlIGZvciB4IGF4aXMgbGFiZWxzKQogIHRoZW1lKHN0cmlwLnRleHQgPSBlbGVtZW50X2JsYW5rKCksIHN0cmlwLmJhY2tncm91bmQgPSBlbGVtZW50X2JsYW5rKCksCiAgICAgICAgcGxvdC5tYXJnaW4gPSB1bml0KGMoMSwgMSwgMiwgMS41KSwgImNtIiksIHBhbmVsLnNwYWNpbmcgPSB1bml0KDEsICJjbSIpKQpgYGAKCmBgYHtyLCBmaWcud2lkdGg9NiwgZmlnLmhlaWdodD04LCBjYWNoZT1UUlVFfQojIHJpZ2h0LXNpZGUgdml6IGRlcGVuZHMgb24gdGhlIG51bWJlciBvZiBzaXRlcyBwZXIgc3BlY2llczoKIyAxIHNpdGUgPSB2ZXJ0aWNhbCBjcm9zc2JhciBvbmx5CiMgMisgc2l0ZXMgPSBwb2ludHMgKyBjcm9zc2JhciBhdCBtZWRpYW4KIyBYKyBzaXRlcyA9IGRlbnNpdGllcyAoY3VycmVudGx5LCBYID0gNCBqdXN0IHRvIGlsbHVzdHJhdGUpCgojIGRpcnR5IGhhY2s6IHggaW4gZnJvbnQgb2YgIlNhbXBsZSBzaXplIiBpcyB0byBoYXZlIHRoYXQgcGFuZWwgc29ydCB0byB0aGUgcmlnaHQgKGFscGhhYmV0aWNhbGx5KSB1bnRpbCBJIGZpZ3VyZSBvdXQgd2h5IGl0IGRvZXNuJ3QganVzdCBnbyBieSBvcmRlci4gVGhpcyBjcm9wcGVkIHVwIGFzIGFuIGlzc3VlIHdoZW4gSSBhZGRlZCB0aGUgZHVtbXkgcG9pbnQgZm9yIHRoZSB4LWF4aXMgZXhwYW5zaW9uLi4uCgojIEFERCBSSUdIVCBGQUNFVApxICU+JSAKICAjIGRlbnNpdGllcyBmb3Igc3BlY2llcyB3aXRoIGVub3VnaCBzaXRlcwogIGZhY2V0X3Bsb3QoInhTYW1wbGUgc2l6ZSIsIGQzYSwgZ2VvbV9kZW5zaXR5X3JpZGdlcywgCiAgICAgICAgICAgICBhZXMoeCA9IG51bSwgZ3JvdXAgPSBsYWJlbCwgZmlsbCA9IGdyb3VwLCBoZWlnaHQgPSAuLmRlbnNpdHkuLiksCiAgICAgICAgICAgICBhbHBoYSA9IC41LCBsd2QgPSAuMywgc2NhbGUgPSAuMykgJT4lCiAgIyB2ZXJ0aWNhbCBjcm9zc2JhciBmb3IgTWRuLCBYIGZvciBNZG4gb2Ygc2l0ZSBtZWRpYW5zCiAgZmFjZXRfcGxvdCgieFNhbXBsZSBzaXplIiwgZDQsIGdlb21fY3Jvc3NiYXJoLCBhZXMoeCA9IE1kbiwgeG1pbiA9IE1kbiwgeG1heCA9IE1kbiwgZ3JvdXAgPSBsYWJlbCwKICAgICAgICAgICAgIGNvbCA9IGdyb3VwKSwgYWxwaGEgPSAuNSwgd2lkdGggPSAuNiwgZmF0dGVuID0gMS41KSAlPiUKICBmYWNldF9wbG90KCJ4U2FtcGxlIHNpemUiLCBkNCwgZ2VvbV9wb2ludCwgYWVzKHggPSBNZG5fc2l0ZTIsIGdyb3VwID0gbGFiZWwpLCBzaGFwZSA9IDQsIHNpemUgPSAyLjUsCiAgICAgICAgICAgICBzdHJva2UgPSAxLjA1KSAlPiUKICAjIHZlcnRpY2FsIG1hcmsgZm9yIGluZGl2aWR1YWwgc2l0ZXMKICBmYWNldF9wbG90KCJ4U2FtcGxlIHNpemUiLCBkM2IsIGdlb21faml0dGVyLCBhZXMoeCA9IG51bTIsIGdyb3VwID0gbGFiZWwpLCBzaGFwZSA9ICJ8Iiwgc2l6ZSA9IDIuNSwKICAgICAgICAgICAgIHdpZHRoID0gLjEsIGhlaWdodCA9IDAsIGFscGhhID0gLjgpCmBgYAoKYGBge3J9Cmdnc2F2ZSgiLi4vZ3JhcGhzL3BoeWxvX3JpZGdlX3NpdGUucGRmIiwgd2lkdGggPSA2LCBoZWlnaHQgPSA4LCBzY2FsZSA9IDIpCmBgYAoKYGBge3IsIGZpZy53aWR0aD02LCBmaWcuaGVpZ2h0PTgsIGNhY2hlPVRSVUV9CiMgQUREIFJJR0hUIEZBQ0VUCnEgJT4lIAogIGZhY2V0X3Bsb3QoInhTYW1wbGUgc2l6ZSIsIGZpbHRlcihkMywgbnVtIDwgMjAwKSwgZ2VvbV9ib3hwbG90aCwgYWVzKHggPSBudW0sIGdyb3VwID0gbGFiZWwsIAogICAgICAgICAgICAgZmlsbCA9IGdyb3VwKSwgYWxwaGEgPSAuNSwgbHdkID0gLjMsIHZhcndpZHRoID0gVCkgJT4lIAogICMgdmVydGljYWwgY3Jvc3NiYXIgZm9yIE1kbgogIGZhY2V0X3Bsb3QoInhTYW1wbGUgc2l6ZSIsIGQ0LCBnZW9tX2Nyb3NzYmFyaCwgYWVzKHggPSBNZG4sIHhtaW4gPSBNZG4sIHhtYXggPSBNZG4sIGdyb3VwID0gbGFiZWwsCiAgICAgICAgICAgICBjb2wgPSBncm91cCksIGFscGhhID0gLjUsIHdpZHRoID0gLjYsIGZhdHRlbiA9IDEuNSkKYGBgCgpgYGB7cn0KZ2dzYXZlKCIuLi9ncmFwaHMvcGh5bG9fcmlkZ2VfYm94LnBkZiIsIHdpZHRoID0gNiwgaGVpZ2h0ID0gOCwgc2NhbGUgPSAyKQpgYGAKCgojIFNlc3Npb24gaW5mbwoKYGBge3J9CnNlc3Npb25JbmZvKCkKYGBgCgo=